Auxiliary-field quantum Monte Carlo study of TiO and MnO molecules 
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Calculations of the binding energy of the transition metal oxide molecules TiO and MnO are 
presented, using a recently developed phaseless auxiliary field quantum Monte Carlo approach. 
This method maps the interacting many-body problem onto a linear combination of non-interacting 
problems by a complex Hubbard-Stratonovich transformation, and controls the phase/sign problem 
with a phaseless approximation relying on a trial wave function. It employs random walks in 
Slater determinant space to project the ground state of the system, and allows use of much of the 
same machinery as in standard density functional theory calculations, such as planewave basis and 
non-local pseudopotentials. The calculations used a single Slater determinant trial wave function 
obtained from a density functional calculation, with no further optimization. The calculated binding 
energies are in good agreement with experiment and with recent diffusion Monte Carlo results. 
Together with previous results for sp-bonded systems, the present study indicates that the phaseless 
auxiliary field method is a robust and promising approach for the study of correlation effects in real 
materials. 

PACS numbers: 02.70.Ss, 71.15.-m, 31.25.-v 



I. INTRODUCTION 

A phaseless auxiliary field (AF) quantum Monte Carlo 
(QMC) method was recently introduced [1] to study cor- 
relation effects in real materials, which has yielded results 
for a variety of sp-bonded materials in good agreement 
with experiment and comparable to those obtained using 
the standard diffusion Monte Carlo method (DMC) [2]. 
In this paper we present the first application of the phase- 
less AF QMC method to the more highly-correlated tran- 
sition metal oxide systems. Because of their complexity 
(the presence of both localized and itinerant characters 
in the electronic degrees of freedom, strong electron-ion 
pseudopotentials, and the presence of many highly cor- 
related electrons), there have been relatively few QMC 
calculations of any type for transition metal systems [3- 
6]. 

There are many important applications based on the 
magnetic, ferroelectric, and superconducting properties 
of transition metal oxides. These effects arise from the 
presence of d-shell electrons whose interactions are of- 
ten highly correlated. The generally successful ah initio 
density functional theory (DFT) approach [7] has had 
limited success in describing these properties, often pre- 
dicting incorrect ground states {e.g. metallic instead of 
insulating). Even in cases where correlation effects are 
less pronounced and the method is qualitatively correct, 
the results are sometimes not of sufficient accuracy. For 
example in ferroelectrics such as PbTiOa , which have es- 
sentially no occupied c?-states, the relatively small and 
usually acceptable DFT errors (~3%) in predicted equi- 
librium volumes can lead to suppression of the ferroelec- 
tric ground state. There is thus a great need for better 
theoretical modeling of transition metal systems. 



Ah initio quantum Monte Carlo methods are an at- 
tractive means to treat explicitly the interacting many 
fermion system. These methods in principle scale alge- 
braically as a low power with system size. However, ex- 
cept for a few special cases, QMC methods are plagued 
by the fermion sign problem [8, 9], which, if uncontrolled, 
results in exponential scaling. No formal solution has 
been found for this problem, but approximate methods 
have been developed that control it. The most estab- 
lished QMC method is the real space fixed-node diffu- 
sion Monte Carlo [10], which has been applied to calcu- 
late many properties of solids and molecules [2]. Recent 
DMC studies have addressed transition metal systems 
such as the TiC molecule [3], TiO and MnO molecules 
[4], solid MnO [5], and solid NiO [6]. 

The new phaseless AF QMC approach [1] is an alter- 
native that has several appealing features. For exam- 
ple, it is formulated in a Hilbert space spanned by some 
fixed one-particle basis, and the freedom to choose any 
one-particle basis suitable for a given problem could be 
advantageous. Moreover, the AF QMC methodology can 
take full advantage of well-established techniques used by 
independent-particle methods with the same basis. With 
a planewave basis, for example, algorithms based on fast 
Fourier transforms (FFT) and separable non-local pseu- 
dopotentials can be carried over from DFT planewave 
codes. Given the remarkable development and success 
of the latter [11], it is clearly desirable to have a QMC 
method that can use exactly the same machinery and 
systematically include correlation effects by simply build- 
ing stochastic ensembles of the independent-particle so- 
lutions. 

The central idea in standard AF QMC methods [12, 13] 
is the mapping of the interacting many-body problem 
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into a linear combination of non-interacting problems in 
external auxiliary fields. Averaging over different AF 
configurations is then performed by Monte Carlo (MC) 
techniques. However, except for special cases (e.g., the 
Hubbard model with on-site interactions), the two-body 
interactions will require auxiliary fields that are complex. 
As a result, the single-particle orbitals become complex, 
and the MC averaging over AF configurations becomes an 
integration over complex variables in many dimensions, 
and a phase problem occurs. 

The phaseless AF QMC method [1] used in this pa- 
per controls the phase/sign problem in an approximate 
manner using a trial wave function. As in fixed-node 
DMC, the calculated results approach the exact ones as 
the trial wave function is improved. The ground-state 
energy in the phaseless method, however, is not a varia- 
tional upper bound [1, 14]. Previous results on sp-bonded 
systems [1, 15] and our current results suggest that the 
calculated energy is quite insensitive to the trial wave 
function. Accurate ground-state energies have been ob- 
tained with simple trial wave functions, namely single 
Slater determinants from DFT or Hartree-Fock calcula- 
tions. 

In this paper, we study the transition metal oxide 
molecules TiO and MnO, using the phaseless AF QMC 
method with planewaves and pseudopotentials. This rep- 
resents the first application of AF-based QMC to tran- 
sition metal oxides. As in regular DFT calculations, 
molecules can be studied with planewaves by placing 
them in large cells (supercells) and using periodic bound- 
ary conditions. This is somewhat disadvantageous be- 
cause one has to insure that the supercells are large 
enough to eliminate the spurious interactions between 
the images of the molecule. Consequently the compu- 
tational cost for isolated atoms and molecules is higher 
than with a localized basis. However, the main motiva- 
tion of the present study is to test the phaseless AF QMC 
method for strongly correlated systems such as transition 
metal oxides, using the same methodology as previously 
used for sp-bonded materials. In addition, a converged 
planewave basis, which is straightforward to achieve aside 
from the computational cost, gives an unbiased represen- 
tation of the Hamiltonian, and facilitates direct compar- 
ison with experiment. 

The remainder of the paper is organized as follows. 
The phaseless AF QMC method is briefiy reviewed in 
Sec. II. The specific formulation using a single-particle 
planewave basis with non-local pseudopotentials is then 
discussed in Sec. HI. Finally, in Sec. IV we present results 
of our calculations for the binding energies of TiO and 
MnO, which are in excellent agreement with experiment. 



II. FORMALISM 

The Hamiltonian for a many-fermion system with two- 
body interactions can be written in any one-particle basis 
in the general form 

M I ^ 

H = Hi+H2 = J2^ijcUj + 2 E Vijki4c]ckci, (1) 

where M is the size of the chosen one-particle basis, and 
cj and Ci are the corresponding creation and annihilation 
operators. Both the one-body (T^) and two-body matrix 
elements {Vijki) are known. 

As in other QMC methods, the auxiliary field quan- 
tum Monte Carlo obtains the ground state \'^g) of H 
by projecting from a trial wave function I^I/t), using the 
imaginary-time propagator e~'^^: 

|*g)oc \im{e-^"r\^T). (2) 

n—*OG 

The trial wave function \'^t}, which should have a non- 
zero overlap with the exact ground state, is assumed to 
be in the form of a single Slater determinant or a linear 
combination of Slater determinants. 

Using a second-order Trotter breakup, we write the 
propagator as: 

The two-body part of the propagator can be written 
as an integral of one-body operators by a Hubbard- 
Stratonovich transformation [16]: 

-'* = Il(;^/>-'"-^-'^'".*<'| w 

after H2 is turned into a sum of squares of one-body 
operators: H2 = —^J2a -^a^a' with Aq, a real number. 

The propagator of Eq. (3) can now be expressed in the 
general form: 

e-''" = j P{a)B{a) da, (5) 

where we have introduced the vector representations a = 
{cti, (72, • • •}, P{<j) is the normal distribution with mean 
and standard deviation 1, and 

B{a) = e-^"^'^ e^"-"- e-^"''^, (6) 

with V = {\/Ai \/A2 W2, • • •}• 

Monte Carlo methods can be used to evaluate the 
multi-dimensional integral of Eq. (5) efficiently. We fol- 
low the procedure [1, 17, 18] of turning the MC process 
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into an open-ended random walk (instead of Metropolis 
sampling of entire paths along imaginary time [12, 13]), 
because it facilitates the imposition of local constraints 
to deal with the sign/phase problem [18]. Each step in 
the random walk takes a Slater determinant to a new 
determinant j^'): 



\<f>'ia))=Biam, 



(7) 



where a is sampled from P{(t). Given sufficient propaga- 
tion time one obtains a MC representation of the ground 

state: |*G> =E<^l</'>- 

This straightforward approach, however, will generally 
lead to an exponential increase in the statistical fluctua- 
tions with the propagation time. One can easily under- 
stand the source of this by realizing that the one-body 
operators v are generally complex, since Aq, usually can- 
not all be made positive in Eq. (4). As a result, the 
orbitals in will become complex as the propagation 
proceeds. This is the phase problem referred to earlier. 
It is of the same origin as the sign problem that occurs 
when B{a) is real. The phase problem is more severe, 
however, because for each \<p), instead of a +\<l)) and 
— symmetry [17], there is now an infinite set {e*^|0)} 
{9 G [0, 2tt)) from which the MC sampling cannot dis- 
tinguish. At large propagation time, the phase of each 
becomes random, and the MC representation of |^'g) 
becomes dominated by noise. 

In Ref. [1] the phaseless auxiliary field QMC method 
was presented to control the phase problem. The first in- 
gredient of this method is an importance-sampling trans- 
formation using a complex importance function, (\E't|^), 
where \^t) is a trial wave function. In the resulting ran- 
dom walk, a walker is propagated to a new position 
\<p') in each step by 



\<f>'ia))=B{a-am. 



(8) 



As in Eq. (7), a is sampled from P{a), but the propagator 
is modified to include a force bias, or shift [19]: 



(1't|v|'^) 
(*t|<^) 



(9) 



A walker carries a weight which is updated according 
to 



w^' = W{(i>) 



(10) 



where W{<i}) can be expressed in terms of the so-called 
local energy, iJ^: 



W^(^) = exp 



exp[-T^i(.^)]. (11) 



In the limit of an exact \^t)i El is a real constant, 
the weight of each walker remains real, and the mixed 



estimate for the energy is phaseless: 



Er- 



(*t|«'g) 



(12) 



With a general \^t) which is not exact, a natural ap- 
proximation is to replace El in Eq.'s (11) and (12) by 
its real part, ReE'i, leading to a phaseless formalism for 
the random walk, with real and positive weights. 

The second ingredient in the phaseless method involves 
a projection: the modified random walk is still "rotation- 
ally invariant" in the complex plane defined by {^t\4>)- 
With the propagation, the walkers will populate the com- 
plex plane symmetrically independent of their initial po- 
sitions. In particular, a finite density of walkers will de- 
velop at the origin where the local energy -E'l(^) diverges, 
and this causes diverging fluctuations in the weights of 
walkers. 

This problem, which is inherent of the "two- 
dimensional" nature of the random walk in the complex 
plane, can be controlled with an additional approxima- 
tion, in which the random walk is projected to "one- 
dimension." This is done, e.g., by multiplying the weight 
of each walker in each step by max{0, cos(A0)}, where 
is the phase of {'^ tW) / t\4') ■ The projection ensures 
that the density of walkers vanish at the origin. Note 
that the projection has no effect when v is real. This 
additional approximation and the importance-sampling 
procedures of Eq.'s (8) through (11) form the basis of 
the new phaseless AF QMC method. 



III. IMPLEMENTATION WITH PLANEWAVES 

The calculations reported in this paper were carried 
out in supercells using a planewave basis and periodic 
boundary conditions (PBC). Pseudopotentials are used 
as in DFT calculations to represent the electron-ion in- 
teraction, eliminating the core electrons from the Hamil- 
tonian. The basis set consists of the M planewaves with 
kinetic energy |k|^/2 < i?cut, where the parameter E'cut 
is a cutoff energy. 

In a planewave basis, the one-body operator Hi of 
Eq. (1) is the sum of the kinetic energy and the electron- 
ion interaction, and H2 represents the electron-electron 
interaction. These can be expressed as: 

2m 



Hi 



k,k',5 

Vjvz,(k,k')4_gCk',s 



k,k',. 



Hn 



E E 



E 

k,k',s,s' ^7^0 



47re^ 



*^k+q,s*^k' 



o/Ck', s'Ck s* 



(13) 
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TABLE I: A summary of the binding energy BE (in eV), equi- 
librium bond length Re (in a.u.) and harmonic vibrational 
frequency co (in cm~^) of the TiO molecule with two different 
pseudopotentials. The first, with Scut = 50 Ry (50 Ry psp), 
was used in all ensuing DFT and QMC calculations. The 
second has a 64 Ry cut-off. The corresponding values of the 
cut-off radius, rc, are listed in the footnotes (in units of a.u.). 
DFT results from both Perdew-Burke-Ernzerhof (PBE) GGA 
[27] and Perdew-Wang 92 LDA [28] functionals are shown, to- 
gether with experimental values. 





BE 


Re 




Experiment [20, 21] 


6.87 or 6.98 


3.06 


1009 


50Rypsp'' GGA 


8.00 


3.02 


1005 


LDA 


9.11 


2.99 


1040 


64Rypsp'' GGA 


7.96 


3.04 


1008 


LDA 


9.05 


3.02 


1044 



"O rc. 1.45 (s), 1.55 (p); Ti r^: 1.40 (s), 1.40 (p), 1.80 (d). 
'O rc: 1.30 (s), 1.39 (p); Ti rc: 1.35 (s), 1.35 (p), 1.52 (d). 



TABLE II: A summary of the binding energy BE (in eV), 
equilibrium bond length Rc (in a.u.), and harmonic vibra- 
tional frequency lu (in cm^^) of the MnO molecule with three 
different pseudopotentials. The first, with Ecut = 64 Ry and 
created from the design non-local (DNL) procedure, was used 
in all ensuing DFT and QMC calculations. Two other sets 
are also tested here, with 64 Ry and 82 Ry cut-off values and 
without DNL. The corresponding Vc values (in a.u.) are listed 
in the footnotes. Calculated results are from DFT GGA. 





BE 


Rc 




Experiment [21] 


3.70 


3.11 


832 


64 Ry DNL psp" 


5.11 


3.11 


822 


64 Ry psp*" 


4.90 


3.07 


878 


82 Ry psp" 


4.99 


3.09 


845 



"O rc 
''O rc 
'^O rc 



1.45 (s), 1.55 (p); Mn rc 
1.45 (s), 1.55 (p); Mn rc 
1.05 (s), 1.02 (p); Mn rc 



1.40 (s), 1.40 (p), 1.65 (d). 
1.40 (s), 1.40 (p), 1.65 (d). 
1.25 (s), 1.25 (p), 1.50 (d). 



IV. RESULTS 



Here c]^ ^ and Ck,s are the creation and annihilation op- 
erators of an electron with momentum k and spin s. 
VL(k — k') and VAri(k, k') arc the local and non-local 
parts of the pseudopotential, respectively, fl is the super- 
cell volume, k and k' are planewaves within the cutoff 
radius, and the q-vectors satisfy |k -|- q|^/2 < £^cut- 

A Hubbard-Stratonovich transformation is applied to 
decouple the electron-electron interaction H2 into a linear 
combination of one-body operators. The resulting one- 
body operators consist of density operators of the form 
= Sk.s ''k4-q,s^k,s- The number of auxiliary fields 
is proportional to the number of unique q vectors that 
the basis allows, i.e., roughly eight times the number of 
planewaves in the basis. 

Non-local pseudopotentials can be treated exactly 
within the present AF QMC formalism, and the use of 
separable forms leads to the same speed-up achieved in 
planewave DFT calculations [1]. This is to be compared 
with standard real-space DMC calculations where an ad- 
ditional locality approximation [22] is used for non-local 
pseudopotentials that depends on the overall quality of 
the trial wave function |5't). (In contrast, the fixed-node 
approximation in DMC only depends on the position of 
the nodal surface of I^E't)-) In order to minimize errors 
due to the locality approximation, small pseudopoten- 
tial cut-off radii rc tend to be used. This could result 
in harder pseudopotentials than otherwise required by 
transferability considerations. In the AF QMC, the use 
of non-local pseudopotentials with larger values rc (de- 
termined only by transferability requirements) does not 
pose any additional difficulty. 



In this paper, we apply the phasclcss AF QMC method 
to calculate the binding energies of the transition metal 
oxide molecules TiO and MnO. Norm- conserving pseu- 
dopotentials are used, and the non-local part of the 
pseudopotential Vnl is represented using the separable 
Kleinman-Bylander (KB) form [23]. 

To obtain the trial wave function |^'t) for each QMC 
calculation, a DFT calculation with the generalized gra- 
dient approximation (GGA) is carried out with the 
ABINIT [24] program, using the same pseudopotentials 
and planewave basis, j^'r) is then taken as the sin- 
gle Slater determinant formed from the occupied single- 
particle orbitals obtained from this DFT calculation, 
with no further optimization. The random walkers are 
all initialized to j^'r), so the many-body ground-state 
projection initiates from the GGA state. In addition, 
this I^'t) is used in the QMC calculations to control the 
sign/phase problem as described in Section II. 

The pseudopotentials were generated by the OPIUM 
program [25] using Ti++, Mn++, and neutral oxygen as 
reference configurations. The titanium and manganese 
semicore states (3s^3p®) were included as valence states, 
so the Ti and Mn atoms contribute 12 and 15 valence 
electrons, respectively, while the O atom contributes 6 
electrons. 

Well-converged planewave cutoffs were 50 Ry for oxy- 
gen and titanium, and 64 Ry for manganese. These Scut's 
were chosen such that the resulting cut-off errors, sys- 
tematically analyzed using DFT calculations, were much 
smaller than the expected QMC statistical errors. In ad- 
dition, we have carried out QMC calculations on a 1 x 1 x 1 
TiO solid supercell with a 50 Ry and 60 Ry cutoff, re- 
spectively. The calculated energies are the same within 
statistical error bars (« 0.1 eV), confirming basis conver- 
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TABLE III: A comparison between LAPW and pscudopoton- 
tial calculations for non spin-polarized TiO in a 7x 7x 14 a.u."^ 
supercell. We show the equilibrium bond length Re (in a.u.) 
and harmonic vibrational frequency tj (in cm~^) from DFT, 
using both GGA and LDA. The two OPIUM pseudopotentials 
are the same as those in Table I. 







Re 


CO 


LAPW 


GGA 


3.01 


1057 




LDA 


2.97 


1097 


50 Ry psp 


GGA 


2.96 


1060 




LDA 


2.94 


1095 


64 Ry psp 


GGA 


2.99 


1058 




LDA 


2.97 


il)91 



gence at the correlated level. The Mn pseudopotential is 
created rising the design non-local pseudopotential pro- 
cedure [26]. This enhances the pseudopotential trans- 
ferability by exploiting the fic;xibility contained in the 
separable KB form of the nonlocal pseudopotential. 

The accuracy of the pseudopotentials was examined 
with DFT calculations of binding energies, as well as the 
equilibrium bond length and harmonic vibrational fre- 
quencies. In Tables I and II, we summarize our calcula- 
tions of these properties for different OPIUM pseudopo- 
tentials. In both cases, increasing the hardness of our 
pseudopotentials did not lead to significant changes in the 
calculated properties. We have also done some of these 
calculations using TrouUier-Martins [29] pseudopoten- 
tials with the same cutofi^ radii, and little difference was 
found. Moreover, our LDA results for the bondlengths 
for TiO and MnO, = 2.99 a.u. and 3.05 a.u., are in 
reasonable agreement with the all-electron LDA values 
{Re = 3.020 a.u. and 3.032 a.u.) [30] and those obtained 
with the Hartwigsen-Goedecker-Hutter pseudopotentials 
[30] . The TiO results of Re and w also compare favorably 
with the calculations of Ref. [31]. 

As a further check on the pseudopotentials, we have 
carried out a comparison between pseudopotential and 
all-electron LAPW calculations. The latter is compu- 
tationally more costly, so we limited the comparison to 
a 7 X 7 X 14 a.u.^ supercell for non spin-polarized TiO 
molecule. Our results for the calculated equilibrium bond 
length and angular frequency of vibration are summa- 
rized in Table III. The close agreement between the 
LAPW and the pseudopotential results gives further ev- 
idence on the reliability of the pseudopotentials. 

Clearly these tests on the quality of the pseudopoten- 
tials are far from perfect. Our pseudopotentials are all 
DFT-bascd, and the tests are with DFT calculations. For 
sj?-bonded systems, we have done plane- wave Hartree- 
Fock (HF) calculations using OPIUM DFT pseudopo- 
tentials, and compared with all-electron HF results. In 
general, these tend to be quite consistent with the DFT 
tests, and often good agreement at the HF level is found 



TABLE IV: A summary of the calculated binding energy of 
the molecule TiO for different supercells. Supercell dimen- 
sions are given in a.u. and binding energies are in eV. The 
QMC statistical errors axe in the last two digits, and axe in- 
dicated in parentheses. At the DFT GGA level, the binding 
energy converges to 8.00 eV. 





GGA 


QMC 


10x11x17 


7.46 


6.59(20) 


12x12x15 


7.77 


6.98(21) 


14x14x15 


7.94 


7.08(21) 


oo 


8.00 





TABLE V: The calculated total ground-state energy of Mn 
for different supercells. Supercell dimensions are in a.u. and 
energies are in eV. The QMC statistical errors axe in the last 
digit, and axe shown in parentheses. 





GGA 


QMC 


11x12x15 


-2766.66 


-2766.40(5) 


12.55x13.69x17.11 


-2766.38 


-2765.66(4) 


14x14x15 


-2766.32 


-2765.89(9) 


15.4x15.4x16.5 


-2766.25 


-2765.74(8) 


oo 


-2766.20 





when good test results have been obtained from DFT 
calculations. Of course, the suitability of a DFT or HF 
pseudopotential (i.e., derived from independent-particle 
procedures) for many-body calculations is a separate is- 
sue, which our tests do not address. Empirically, such 
pseudopotentials have been widely used in many-body 
calculations and have been quite successful. 

The use of PBC with a planewave basis requires su- 
percells that are large enough to control spurious inter- 
actions between the periodic images of the system under 
study. We studied convergence with respect to such size 
effects using both ABINIT and QMC calculations. Rep- 
resentative results are shown in Tables IV and V. 

Estimating the size-effects in the AF QMC calculations 
is complicated by the presence of finite Totter time-step 
(r) errors. The QMC values shown in Tables IV and V 
are final values after extrapolations in r, as discussed fur- 
ther below. The range of supercells shown in Table IV 
corresponds to about 12,000-17,000 planewaves in our ba- 
sis. For the Ti atom, the largest two supercells resulted in 
a degeneracy of the highest-lying occupied d-orbitals in 
the density functional calculations. To break the degen- 
eracy, these supercells were modified to 11.6 x 12 x 15 a.u.'^ 
and 13.5 x 14x 15 a.u.^, respectively. The fully converged 
value of the DFT GGA TiO binding energy is 8.00 eV, 
as shown. For the AF QMC calculations, the binding 
energies for the larger sizes are converged to well within 
the statistical errors. 

Table V shows the energy of the Mn atom for different 
supercell sizes. The corresponding number of planewaves 
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-2007.5 
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<a 



2008.0 



-2008.5 



-2009.0 




FIG. 1: QMC time-step r dependence of the total energy of 
the TiO molecule. A 10 x 11 x 17 a.u.^ supercell was used. The 
solid line is a linear fit to the calculated QMC energies (solid 
squares). The final extrapolated energy E = —2007.72(17) eV 
is shown as a star. 



is between 17,000 and 34,000. As can be seen, the QMC 
energy is converged to less than the statistical error for 
the 14 X 14 X 15 supercell, ahhough for the smaller su- 
percells, the finite-size errors are significant both in GGA 
and QMC. The MnO molecule, on the other hand, ex- 
hibits a much smaller size effect, with QMC energies of 
-3195.50(11) eV and -3195.58(7) eV for the 11 x 12 x 15 
and 14 X 14 X 15 supercells, respectively. 

The QMC Trotter errors were examined by studying 
the individual time-step dependence for the atoms and 
the molecule using a particular supercell size. Figure 1, 
for example, shows the Trotter extrapolation for the TiO 
molecule, done with a 10 x 11 x 17 supercell. The Trotter 
behavior obtained from this procedure was then used to 
extrapolate the QMC data of other supercell sizes, for 
which calculations were performed with the time step 
fixed at r = 0.025 Ry~^. The final extrapolated re- 
sults are what is shown (e.g., in Table IV). Figure 2 
shows the time-step dependence of MnO, which exhibits 
a quadratic behavior compared to the more linear depen- 
dence in Fig. 1 for TiO. The Mn and O atoms exhibit 
much smaller finite-r errors, as is also the case with the 
Ti atom (not shown). 

Table VI summarizes the results for the molecular 
binding energies. For comparison we also include results 
from a recent diffusion Monte Carlo study by Wagner 
and Mitas [4]. As mentioned, our AF QMC calculations 
use a single-determinant trial wave functions obtained 
from a DFT GGA calculation, without a Jastrow factor 
or any further optimization to the determinant. We see 
that the calculated binding energies from AF QMC and 
those from DMC [4] with trial wave functions containing 
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-3196 
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-3197 
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FIG. 2: QMC time-step r dependence of MnO, Mn, and O. 
An 11 x 12 X 15a.u.^ supercell was used for MnO and Mn, 
and a 10 X 11 X 17a.u.^ supercell for oxygen. The solid line is 
a least-squares fit to the QMC energies (solid squares). The 
final extrapolated values are shown as a star. MC statistical 
error bars are indicated. 



either an optimized hybrid B3LYP determinant or mul- 
tiple determinants from MCSCF are in good agreement 
with each other and with experiment. DMC with a trial 
wave function containing the Jastrow and a single Slater 
determinant from HF, on the other hand, gives somewhat 
worse agreement with experiment. We have not carried 
out AF QMC calculations using an HF trial wave func- 
tion for these molecules. In several sp-bonded molecules, 
DFT and HF-generated trial wave functions showed little 
difference in the calculated energies in AF QMC. 

We have also included in Table VI the results for the 
binding energy of the O2 molecule. Because of the short 
bond length of this molecule (Re=2.281 a.u.), a harder 
pseudopotential was used, with a higher Ecut of 82 Ry 
and smaller values for Tc (last entry in Table II). At 
the DFT GGA level the binding energy is 5.72 eV. Our 
QMC results shown in Table VI were obtained using a 
supercell of size 8 x 9 x 11 a.u.*^. Additional QMC and 
DFT calculations with a larger supercell of 11 x 12 x 
13 a.u.'^ have verified that the finite-size effects are within 
our statistical error bars (« 0.1 eV). Again, we see that 
the agreement with experiment is very good. 
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TABLE VI: A summary of the binding energies of the 
molecules TiO, MnO, and O2. Calculated results from the 
present QMC method and diffusion Monte Carlo (TiO and 
MnO from Ref. [4], and O2 from Ref. [32]) are shown, together 
with experimental values (TiO from Refs. [20, 21], MnO from 
Ref. [20], and O2 from Ref. [32]). Equilibrium experimen- 
tal bond lengths were used in the molecule calculations. Our 
QMC used as trial wave function a single Slater determinant 
from DFT GGA. The trial wave functions used in DMC are 
indicated in the footnotes. All energies are in eV, and the 
experimental zero point energy is added to each molecule. 





TiO 


MnO 


O2 


Experiment 


6.98 


3.70 


5.1152(9) 




6.87 






DMC (HE)" 


6.3(1) 


2.9(1) 


4.84(2) 


DMC (B3LYP)'' 


6.9(1) 


3.4(2) 




DMC (MCSCF)<= 


6.7(2) 


3.4(2) 




Present QMC 


7.02(21) 3.79(34) 5.12(10) 



"Trial w.f. is a (HF 1-dctcrminant) X Jastrow. 

'Trial w.f. is a (DFT B3LYP l-determinant)x Jastrow. 

"^Trial w.f. is a (MCSCF multi-determinant) x Jastrow. 



putationally rather demanding, especially with transition 
metal oxides. For instance, the ground-state energy of 
the MnO molecule in Fig. 2 at the single Trotter step 
of T = 0.008 Ry^^ (with an error bar of 0.35 eV) was ob- 
tained from running on an Intel XEON cluster (3.2 GHz) 
for about 150 hours using 72 processors. 

In summary, we have presented the first study of tran- 
sition metal oxide molecules by AF QMC. Wc have shown 
that the binding energies of TiO and MnO calculated 
with the new phaseless AF QMC method [1] are in good 
agreement with experiment, and arc comparable to the 
best results obtained from diffusion Monte Carlo [4]. It is 
encouraging that a trial wave function of only DFT single 
Slater determinants was sufficient for the phaseless QMC 
method to reach this accuracy. Together with previous 
results for sp-bonded systems [1, 15], the present study 
indicates that the phaseless method is a robust QMC 
method. Complementary to standard DMC, it offers a 
promising approach for the computation of correlation 
effects in real materials. 

V. ACKNOWLEDGMENTS: 



Finally, we comment briefly on the computational 
cost. As mentioned, the use of planewaves for iso- 
lated molecules is somewhat disadvantageous even at the 
density-functional level, because of the need for large su- 
perccll sizes to reduce the spurious interactions between 
the images of the molecule. The number of planewaves, 
M, is proportional to the supercell volume, and the com- 
putational cost scales with M as MlnM. (In addition, 
it scales quadratically with the number of electrons.) As 
a result, these planewave AF QMC calculations are com- 
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